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Abstract. In the theory of diffusive acceleration at 
oblique shock fronts the question of the existence of a dis- 
continuity of energetic particle density is contentious. The 
resolution of this problem is interesting from a theoreti- 
cal point of view, and potentially for the interpretation 
of observations of particle densities at heliospheric shocks 
and of high-resolution radio observations of the rims of 
supernova remnants. It can be shown analytically that an 
isotropic particle distribution at a shock front implies con- 
tinuity of the particle density - whether or not the shock 
is oblique. However, if the obliquity of the shock induces 
an anisotropy, a jump is permitted. Both semi-analytic 
computations and Monte-Carlo simulations are used to 
show that, for interesting parameter ranges, a jump is in- 
deed produced, with accelerated particles concentrated in 
a precursor ahead of the shock front. 

Key words: Acceleration of particles - Shock waves - 
Methods: numerical - cosmic rays - ISM: supernova rem- 
nants 



1. Controversy — general description 

In the test-particle theory of diffusive shock acceleration, 
the phase-space spectral index s of accelerated particles 
depends solely on the compression ratio r — p' j p of the 
shock (where p and p' are the upstream and downstream 
densities respectively): s = Zr/{r — 1), which results in 
s = 4 for a strong shock in an ideal gas with Cp/cv — 5/3 
(Axford et al. |1977| ; Krymskii |l977t Bell |1978| ; Blandford 
& Ostriker 197^ ). This result, like many other analyti- 
cal predictions (for a review see Drury 1983; Blandford & 
Eichler 1987), depends on the assumption that the phase 



space density is close to being isotropic, even at the shock 
front. In this case, the density profile of accelerated test 
particles is a continuous function of position (e.g.. Kirk 
et al. 1994, page 262). In planar symmetry, the density 
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is constant downstream and drops off exponentially up- 
stream of the shock. This situation applies to both parallel 
and obli que s hocks, as was pointed out by, for example, 
Axford ( |l98lD . 

Recently, a discussion has arisen in the literature con- 
cerning the occurrence of discontinuities in the density of 
accelerated part icles at an oblique shock front. Whereas 
Ostrowski ( 199l|) (hereaft er 09 1) finds a substantial ef- 
fect, Naito & Takahara ( |1995| ) (hereafter N&T95) as- 
sert that the density is continuous. Both of these papers 
present Monte-Carlo simulations of particle acceleration 
in which the velocity of the shock front is a significant 
fraction of the speed of light, and take explicit account of 
a possible anisotropy of the particle distribution. From the 
theory of diffusion, it is well-known that the anisotropy of 
particles of speed v is of the order of u/v, where u is the 
speed at which the shock sweeps through the medium re- 
sponsible for making the particles diffuse. At an oblique 
shock front, the speed relevant for particles diffusing along 
a magnetic field line is the velocity of the intersection point 
of the shock front and a given field line, u = Wg/cos*!', 
where Ug is the inflow speed along the shock normal and 
$ is the angle between the shock normal and the magnetic 
field, measured in the upstream rest frame of the plasma. 
Thus, if the field is oblique, even a relatively slow shock 
front can produce a substantial anisotropy. The question 
which we address in this paper is whether or not this 
anisotropy leads to a discontinuity in the particle den- 
sity and under what conditions such an effect could be 
observed. 

An approximation which is often used in treating 
oblique shocks is that in which a particle crossing the 
shock conserves its magnetic moment (also referred to as 
the first adiabatic invariant; e.g., Webb et al. 1983| ). This 
approximation is valid for non-relativistic perpendicular 
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shocks (Whipple et al. 1986). Terasawa ( |1979 ) and Decker 
( |1988D have performed numerical simulations which con- 
firm its accuracy at non-relativistic oblique shocks and 
Begelman & Kirk ( |1990| ) have conducted tests at rela- 
tivistic shocks. It is important to bear in mind, however, 
that this approximation breaks do wn for sufficiently fast 
shocks (see, for example. Kirk et al. 1994 , page 241). Con- 
servation of magnetic moment implies that a particle can 
be reflected by the magnetic compression at a fast-mode 
shock front, and the question of the existence of a discon- 
tinuity in the particle density is intimately connected with 
the phenomenon of reflection, as pointed out in an early 
paper on this subject (Achterberg & Norman 1980| ). If this 
approximation is adopted, the problem of acceleration at 
an oblique shock of particles which undergo pitch-angle 
diffusion along field lines can be solved semi-analytically 
(Kirk & Heavens |1989| , hereafter K&H89), at least for the 
case in which the accelerated particles are ultra-relativistic 
(v = c) and build a power-law distribution in momentum. 
In Sect. 0, we address the question of the existence of a 
discontinuity using a new Monte-Carlo code incorporat- 
ing the assumption of conservation of magnetic moment 
(Gieseler 1998). Since there already exist contradictory 
simulations in the literature, we take particular care to 
check the results obtained with this new code against the 
analytic results of K&H89. Both the code and the anal- 
ysis find a discontinuity. However, although N&T95 also 
adopted the assumption of conservation of magnetic mo- 
ment in their simulations, they did not find a discontinuity. 
We show that this can be caused by low spatial resolu- 
tion of the simulation and does not necessarily invalidate 
results con cernin g the spectrum of accelerated particles. 
Ostrowski (1991), on the other hand, did find a disconti- 
nuity, without implementing the conservation of magnetic 
moment. In Sect, g we argue that Liouville's theorem in 
fact guarantees continuity in this case. This result can once 
again be caused by a rapidly varying density profile which 
is sampled at low spatial resolution; in this case it appears 
because the length scale of a gyration radius was not re- 
solved. However, Ostrowski's results do indicate that the 
accumulation of particles upstream of a relativistic shock 
front is an important physical effect, even though there 
may not be a discontinuity in the strict sense. 

At non-relativistic shocks, such as those observed in 
the solar system, the analytic method of K&H89 fails, 
and one must rely on simulations. Here again the situ- 
ation is not clear-cut. Monte-Carlo simulations by Ellison 
et al. ( 1996D do not show discontinuities or a jump across 
the shock front, whereas recent numerical solutions using 
the finite difference met hod have found such effects (Ruf- 
folo 1999). In Sect. 3^ we present high resolution test- 
particle simulations using parameters appropriate for solar 
system shocks and demonstrate the existence of accumu- 
lated particles in the precursor of the shock. 

The resolution of this question is not only of formal 
interest, but is also relevant for the interpretation of data 



taken by the Ulysses spacecraft and of observations of the 
radio emission of supernova remnants. These applications 
are discussed briefly in Sect. H. 

2. Analytical considerations 

The system we consider consists of energetic particles 
which move in two half-spaces separated by a shock front. 
Plasma flows into and out of the shock front carrying with 
it a uniform magnetic field. In addition to the uniform 
field, we assume there exist magnetic fluctuations static 
in the rest frame of the plasma on each side of the shock 
whose effect on the particle motion may be described via a 
scattering operator. Two levels of approximation are im- 
portant. In the first, the particle is described by its full 
six-dimensional phase space coordinates. The trajectory 
is integrated in an explicit realisation of the stochastic 
magnetic field, taking full account of the gyration phase. 
Upon crossing the shock front, the trajectory experiences 
no forces other than those exerted by the magnetic field. 
This is the approach used by 091. The second level of ap- 
proximation is one in which the particle is described by five 
coordinates: the position of the guiding centre, the magni- 
tude of the momentum and the pitch angle. The magnetic 
field is assumed uniform and the guiding centre follows a 
field line except when crossing the shock. The magnetic 
moment is conserved in between scattering events, which 
are assumed to change only the pitch angle, and also upon 
crossing the shock front. This is the approximation used 
by K&H89 and N&T95. It is generally referred to in the 
literature as the 'drift approximation', and we adopt this 
terminology here, although we do not explicitly consider 
the drifts themselves. 

In both levels of approximation a discontinuity can in 
principle arise in the formal description of the dependence 
of the distribution function on the spatial coordinate along 
the shock normal. However, the physical interpretation of 
a mathematical discontinuity is different in each case. 

In the first approximation, the shock front is taken to 
have zero thickness. In reality, the shock transition will 
extend over a finite region in space, which is assumed to 
be small compared to the other length scales of interest, 
for example the gyration radius of the energetic particles. 
Thus, if the acceleration process were to produce a signifi- 
cant change in the energetic particle density over a length 
scale smaller than the gyration radius, but perhaps com- 
parable to the shock thickness, this would appear as a 
discontinuity in simulations such as those of 091. How- 
ever, since particles do not suffer impulsive deflection at 
the shock front, i.e., the momentum coordinates of a par- 
ticle are in general continuous functions of position even 
at the shock front itself, it follows from Liouville's theo- 
rem that the phase-space distribution function, which is 
constant along trajectories, is also a continuous function 
of position across the shock front. The particle density is 
simply an integral of this function over all momenta, so 
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that it is also continuous, provided the momentum is mea- 
sured in the same frame of reference both upstream and 
downstream of the shock. Thus, there can be no formal 
discontinuity in this approach. Of course, the density may 
vary smoothly on the length scale of the gyration radius, 
so that we can interpret Ostrowski's result as due to a 
strong gradient in the density on this length scale. 

The drift approximation, however, can only resolve 
changes which occur on length scales longer than a gy- 
ration radius, so that the distribution found by Ostrowski 
must appear as a discontinuity in simulations which use 
this approximation. There is a close relationship between 
such a discontinuity and the angular distribution, which 
can be understood as follows. Consider an oblique shock 
viewed in the de Hoffmann- Teller (hereafter dHT) frame 
(de Hoffmann & Teller |1950| Kirk et al. |1994D , in which 
the electric field vanishes and the shock is stationary. A 
particle trajectory is now described by only five coordi- 
nates, of which the magnetic moment is conserved both 
between scatterings and on encountering the shock front. 
Denoting by pj^ the component of the particle momentum 
perpendicular to the magnetic field B, the magnetic mo- 
ment is p\/B. Because p = \p\ is also conserved in the 
dHT frame (in which the electric field vanishes) this leads 
to 



1 - /.2 _ 1 _ (^' 



B 



B' 



(1) 



where a prime denotes downstream quantities. Rearrang- 
ing, the downstream pitch angle is given by 



M' = (/./|/i|)v/(M2-A^L)/(l-/^L): 



(2) 



where the cosine of the of the loss-cone angle is given by 
/i„it = a/1 — B/B'. Thus, according to its pitch angle, an 
upstream particle may be reflected or transmitted, and we 
can divide phase space on the upstream side of the shock 
into four regions: 

(a) /icrit < /^ < 1, particles approaching the shock which 
will be transmitted (i.e., in the loss cone); 

(b) < ^ < /icrit, particles approaching the shock which 
will be reflected; 

(c) — /x„it < fj, < Q, particles leaving the shock after reflec- 
tion; 

(d) — 1 < /I < — ^„it, particles leaving the shock after 
transmission from downstream. 

The phase space downstream splits into just two regions, 
/x' > and /x' < 0, since no trajectories incident from 
downstream are reflected. 

In the drift approximation, the phase-space distribu- 
tion function is indep enden t of gyration phase to lowest 
order (e.g., Spatschek 199C , page 145), so that application 
of Liouville's theorem yields 

/(P,A*) = /(P, -/") for 1^*1 < Merit (reflection), (3) 

/(PiA') = /'b, m') for ImI > M„it (transmission), (4) 



where conservation of the momentum p in the dHT frame 
is used. 

The upstream and downstream densities are given by 

-Pcrit Merit 1 

n = / dnnt{fJ.)+ / dfin,:{fi)+ dfint{fi), (5) 

-1 — /^crit Merit 



n' = /dAi'nt(/i') 

-1 



(6) 



where we have defined the quantities 
"-t,r(Ai) := f{p,^J.)2^^p'^dp, (7) 

n[{^l') :^ f f'{p,y.')2np'dp=nt{p), (8) 



such that the suffix r refers to particles which are or will 
be reflected (i.e., |/i| < /icrit) and the suffix t to ones which 
are or will be transmitted (|/z| > //crit, or —1 < /i' < 1). 
The second relation in Eq. (8) follows from Eq. (4). 

A continuous density distribution at the shock front 
{n = n') is in general not guaranteed. This can be 
seen by using specific assumptions about the density 
of trans mitted and reflected particles and Eqs. (5)-(6) 
(Gieseler 199S ). For the physical distribution of Fig. ^, the 
values of the density which can be calculated from Eq. (5) 
and Eq. (6) are quite different, as discussed below. The 
reason for the upstream and downstream densities to be 
different can be understood from Eq. (8), which preserves 
a 'balance' of the transmitted particles, whereas the re- 
flected particles contribute only to the upstream density. 
A continuous density is obtained if there is no compres- 
sion of the magnetic field {B = B'). In this case no par- 
ticles are refiected (/i„it = 0) and Eqs. (5)-(8) then give 
n = n' . This is valid for a parallel shock and the trivial 
case in which no shock front is present. If the pitch-angle 
distribution is isotropic at an oblique shock front, then 
one again finds n — n' hy integration of Eqs. (5) and (6). 
Noting that in this case nt{fi) and n^dj,) have the same 
constant value, continuity follows from Eq. (8). In other 
words, the contribution of the refiected particles upstream 
exactly balances the compression of the transmitted ones 
downstream in this case. 

To summarise this section, in general the density of ac- 
celerated test particles at oblique shocks will vary on the 
length scale of the gyration radius across the shock front. 
This variation is closely connected with the anisotropy of 
the angular distribution. It should appear as a disconti- 
nuity in treatments which use the drift approximation to 
the particle motion. However, when conditions are such 
that the theory of diffusive acceleration applies, i.e., the 
particle velocity is much larger than the shock speed pro- 
jected along the magnetic field {v ^ Ws/cos<I>), then the 
anisotropy and the associated density variation are small. 
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3. Monte-Carlo simulations 

We now present results from test-particle simulations of 
accelerated particles at shock fronts. The key aspects 
of the technique have been used and described by sev- 
era l auth ors (e.g., Kirk & Schneider 1987 ; 091; Baring et 
al. 1993 ; N&T95), so that a brief description suffices. We 
consider oblique shocks, where the magnetic field is in- 
clined at an angle $ with respect to the shock normal 
in the upstream rest frame, and has no dynamical effect 
on the plasma flow. The shock speed in this frame is Ug. 
We consider (in principle) the whole range of <i> for sub- 
luminal shocks, e.g. Ug < Ms/cos<1> < 1 (here and below: 
c = 1). The gyration centre of a particle's trajectory is 
followed in the upstream and downstream rest frames of 
the background plasma. In these frames, as in the dHT 
frame, the momentum p — \p\ is constant. Particles move 
along the magnetic field under the influence of small scale 
irregularities which lead to pitch-angle scattering. We do 
not investigate the effect of transport of particles perpen- 
dicular to the mean magnetic field. We use an algorithm 
for calculating a pitch angle /i„„„ from a given pitch an- 
gle /i which was given by 091 (see Fig. 1 therein). From 
two random numbers i?i and i?2 which are uniformly dis- 
tributed on the interval [0, 1], the new pitch angle is given 

by 

cos An = 1 - (1 - COS Arj,„,J Ri , (9) 

fj,^^^ = ^ cos Ar2 + \/l - /i2 sin An cos(27ri?2) , (10) 

where we have chosen Ail^^^ = 0.1 for most of the simula- 
tions shown below. This gives a very good approximation 
of pitch-angle scattering with an infinitesimal amplitude. 
The results (at least for the spectral index and the density 
distribution) do not change substantially for a factor of 5 
higher or lower value of Afi^^^.E! The time step At for suc- 
cessive scatterings is kept constant. These scatterings are 
performed in the upstream and downstream rest frames. 
At a shock crossing, we transform the particle momentum 
p and the pitch angle /i into the dHT frame, which is al- 
ways possible at subluniinal oblique shocks (for a descrip- 
tion of the Lorentz transformations between these frames 



can be estimated as 



see K&H89 or Gieseler |1998|), and apply the conservation 
of magnetic moment (Eq. 1) to find the new downstream 
pitch angle in the dHT frame. Transformation into the 
downstream rest frame then gives the new values of the 
momentum and pitch angle. This method is also used by 
K&H89 and N&T95, to which we compare our results. 

The validity of the assumption of conservation of the 
magnetic moment can be evaluated as follows. First, con- 
sider a 'strict' condition which can be applied to a single 



particle crossing the shock (Kirk et al. 1994, page 241). 
The number of times a trajectory intersects the shock front 



N.., 



tana tan<l>d 



(11) 



where a = cos""'^ /i is the pitch angle, and $dHT is the an- 
gle between the upstream field and the shock normal mea- 
sured in the dHT frame.a Adiabatic behaviour is expected 
when N^rass ^ 1. An accurate separation between trans- 
mitted and reflected particles is guaranteed if we demand 
this condition to be valid for at least all particles outside 
the loss cone, i.e., \fi\ < fj,„it = cosa^rit = v^l — B/B' , 
where B/B' = [(1 + tan^ $dHT)/(l + r^ tan^ ^^ht)]"^/^ 
and r is the compression ratio. Inserting the critical pitch 
angle adt in Eq. (11) and expanding in powers of l/N^^^^^ , 
we get a condition for the upstream inclination angle in 
the dHT frame which is independent of the fluid speed: 



tan<l>rt 



> Vr 



1. 



(12) 



^ 091 used Af2„ax ~ 0.3 as an approximation for infinitesi- 
mal pitch-angle scattering. 



For the non-relativistic shock speed Us = 0.01, Eq. (12) 
is fulfilled for essentially the whole parameter range in 
Us/ COS $ considered here. For shock speeds in the mildly 
relativistic region, Eq. (12) always holds in an (albeit 
small) region at Us/cos$ — > 1. 

The condition A^croaa ^ 1 is, however, much too strict. 
A less conservative approach is to assume validity of the 
approximation in an average (over phases) sense. This has 
been shown b y simu lations as mentioned above. For ex- 
ample Decker (1988) finds that for non-relativistic shocks 
the results compare well with adiabatic behaviour over 
the whole range of <&. In particular, the data shown in 
his Fig. 8 show close agreement concerning the reflection 
probability for the approximate and numerical results. Os- 
trowski ( 1991 ) used an extended Monte-Carlo code which 
does not assume adiabatic behaviour. He was thus able 
to check the quality of the approximation at mildly rel- 
ativistic shock speeds. His Fig. 3 shows good agreement 
with K&H89 for the spectral index. Departures start to 
occur at Ug = 0.5, but always diminish as Ms/ cos $ — > 1. 
Furthermore, we have directly compared our simulations 
with his for Ug — 0.3 at $ = 60° and 70° and found the 
same density discontinuity (see below). 

To test the accuracy of our scheme against the calcu- 
lation of K&H89, particles are injected in the upstream 
plasma with a velocity corresponding to jv = p/m > 2, 
which means that they are already relativistic at injec- 
tion. For intersection speeds of shock and magnetic field 
of Us/ cos($) > 0.9 the injection momentum has to be 
even higher. We are interested in particles with momenta 
where the distribution has attained a pure power law, 
and especially in the spectral index of this power law. 
Only particles from the power-law region in the momen- 
tum distribution are used for the pitch-angle and density 
distributions. Within this region we found that the an- 
gular distribution is independent of the momentum range 



Note that cos^ $ = (1 + u^ tan^ <l>dHT)/(l -|- tan^ <l>dHT) 
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5.068 ± 0.1818E-02 
0,1 709 ± 0,3590E-03 




18 20 

Log(p/m) 

Fig. 1. Monte-Carlo simulation of a momentum spectrum for 

compression ratio r — 4, shock speed Ug — 0.4 and an angle of 

the magnetic field of <l> = 45°. The speed of the intersection 

of magnetic field and shock is then Us/cos<i> = 0.5657 (same 

parameters as for Fig. | and Fig. ^. The spectral index of the suits from K&H89 (sohd line: 

phase space density is s = 3.1709 ± 0.0004. dotted line; Us = 0.01). 
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Fig. 2. Comparison of the spectral index from Monte-Carlo 
simulations for infiiutesimal pitch-angle scattering with 
Ailmax =0.1 (dots, stars, squares) with semi-analytical re- 
Ms = 0.5; dashed line: Us = 0.3; 



chosen, and the power-law index is thus independent of an- 
gle. Simulation of an individual particle terminates when a 
maximum momentum is reached (two orders of magnitude 
higher than the upper limit of Fig. m), or when the parti- 
cle reaches a distance from the shock on the downstream 
side which is three times greater than the left boundary 
of Fig. H, at which the density distribution has already 
reached its constant downstream value. These limits vary 
with the parameters Ug and $, and were always chosen 
such that boundary effects are not important. The results 
from Monte-Carlo simulations and semi-analytical com- 
putations arc shown (with the exception of Fig. ph for a 
compression ratio of r = 4. 

The three main aspects of the results are spectral in- 
dex, angular distribution and spatial distribution, and 
these are presented in turn. 

3.1. Spectral index 

The momentum distribution of particles is measured at 
the shock front in the upstream rest frame. An exam- 
ple is shown in Fig. |^, where 10^ independent particles 
were simulated. The plot also contains a fit to the function 
y = IQP^+P^'^^ together with the values of the fit param- 
eters and their statistical errors. The spectral index can 
be calculated from AN/Alog{p/m) oc p~'*+'^; this gives 
s — 3 — p2. The statistical error of the fit to the spectra 
is less than 1% for all spectral indices discussed below. To 
achieve maximum accuracy from the Monte-Carlo simula- 
tions, we fit the spectrum over a large finite range of mo- 
mentum and did not include loss mechanisms. For s < 3.2 
the statistical errors are less than 0.1%. For spectral in- 
dices near s ~ 3 the statistics are very much better, be- 
cause particles then gain energy mainly due to reflection. 
This means that not much time is taken in following the 
particle trajectories in those parts of the downstream re- 



gion where the particle has a low chance of returning to the 
shock front. As a result, the error in s is less than 0.01% for 
spectral indices corresponding to Ug/cos'I' > 0.8. Fig. || 
shows the spectral indices for non-relativistic and mildly 
relativistic shock velocities and a wide range of inclination 
angles $. The lines in this plot are taken from Fig. 2 and 
Fig. 6 of K&H89 and are in precise agreement with the 
simulations. 

Very important for an understanding of the spectral 
index and the density profile is the effect of the underly- 
ing scattering law. The results become quite different for 
isotropisation after each scattering, given by Arj^^j^ = tt. 
Although this does not correspond exactly to the sim- 
ulation of large-angle scattering on point-like scattering 
centres (for which it would be necessary to choose expo- 
nentially distributed time steps, and consider cross field 
diffusion), the effects it produces should be qualitatively 
similar. Using a scattering law in the transition regime 
between infinitesimal pitch-angle scattering and isotropi- 
sation after each scattering allows an investigation of this 
dependence. We focus on one example with Ari„,^^ = 1.0 
(intermediate pitch-angle scattering). This kind of scat- 
tering can have a very big influence on the precursor of 
accelerated particles. Under pitch-angle scattering with 
Ari„^x ^ 1 particles have a higher escape probability from 
the shock because they are free to change their pitch angle 
to a value in the loss cone (/Ucrit < A* < 1) within a few scat- 
terings, and may then be transmitted through the shock. 
This reduces the upstream density, as we will see later, 
but also leads to a steeper spectrum for accelerated test 
particles. Figure || shows a comparison of the spectral in- 
dex for intermediate and infinitesimal pitch-angle scatter- 
ing. The stars, squares and dots show Monte-Carlo results 
for intermediate pitch-angle scattering with AfJ^^,,, — 1.0, 
whereas the lines are (as in Fig. g) results from K&H89 for 
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Fig. 3. Comparison of the spectral index for intermediate 
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infinitesimal pitch-angle scattering {Afl^^^ <C 1). It can be 
seen that even for a relatively large value of Afi^^^ = 1.0, 
the effect is very small for most intersection velocities of 
magnetic field and shock. 

For infinitesimal pitch-angle scattering and Ug/ cos $ > 
0.8, a flat s ~ 3 spectrum is achieved. In the standard 
picture, this spectrum corresponds to acceleration with 
vanishing escape probability. Here it is associated with a 
very strong pile-up (see Fig. ^. A large (~ 1) value of 
Afl„,^ reduces the pile-up (see Fig. 0) and effectively in- 
creases the escape probability, leading to the steeper spec- 
trum shown in Fig. for Wg/cos*!' > 0.95. For Ar2,„ax — ti" 
(isotropisation at each scattering) the spectrum becomes 
steeper for all inclination angles. This point was noted 
by N&T95 (their Fig. 7), and our results for the spectral 
index are in good agreement with theirs. 

3.2. Pitch-angle distribution 

As shown in Sect. ||, the pitch- angle distribution of accel- 
erated particles plays a crucial role in determining whether 
the density distribution has a jump at the shock front. 
For every set of parameters we have measured the pitch- 
angle distribution and present an example of the results 
in Fig. H. The lines represent the distributions calculated 
semi-analytically by K&H89, whereas the discrete sym- 
bols show the contents of the 'bins' filled in the above 
described Monte-Carlo method. The four representations 
of the pitch-angle distribution at the shock shown in Fig. H 
are: 1. Upstream of the shock in the upstream rest frame 
(Fig. ^ dashed line and open triangles). 2. Upstream 
of the shock in the dHT frame (Fig. ^ solid line and 
open circles). 3. Downstream of the shock in the down- 
stream rest frame (Fig. ^ dashed line and open triangles 
4. Downstream of the shock in the dHT frame (Fig 
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Fig. 4. Pitch-angle distribution at the shock front for 
Af2,-„i,x =0.1 and compression ratio r = 4. Figs, a) and b) 
show the distribution at the upstream and downstream side of 
the shock. In each figure the lines are from K&H89 and the dis- 
crete symbols display the 'bins' of the Monte-Carlo simulations. 
The solid lines and the circles show the distribution in the dHT 
frame, and the dashed lines and the triangles show the distri- 
butions in the upstream (a) and downstream (b) rest frame. 
See text for the normalisation. The shock speed is Ms = 0.4, and 
the upstream inclination angle of the magnetic field $ = 45° 
(compare Fig. hi and Fig. pt). 



solid line and open circles). The upstream distributions 
(Fig. Ha) are normalised to their maximum value. The 
normalisation at the downstream side of the shock are 
as follows: the Monte-Carlo distribution in the dHT frame 
has the same normalisation as the corresponding upstream 
one, to allow a direct comparison between the two (see be- 
low); however, the Monte-Carlo distribution in the down- 
stream rest frame is normalised to a maximal value of 0.4 
to enable it to be displayed in the same figure. In each 
case, the normalisation of the semi-analytical results of 
K&H89 are chosen to provide the best fit to the Monte- 
Carlo distributions. Comparison of the pitch-angle distri- 
butions from the two methods in Fig. ^ confirms that they 
are in close agreement. 

Particles with pitch angle /i > in the dHT frame 
move in the downstream direction. From Eq. (1) we can 
calculate a critical pitch angle cosine below which these 
particles will be reflected upstream of the shock. For r = 4, 
$ = 45'' and Us = 0.4 we get /x„it = 0-826. Reflected 
particles contribute symmetrically about /x = and con- 
stitute the major part of the distribution (see solid line 
of Fig. |ja), indicating that repeated reflections are ef- 
fective in keeping particles upstream. The particles with 
/i < — /i^Ht = —0.826 in the dHT frame are particles which 
cross the shock from the downstream side. A Lorentz 
transformationtl of the critical angle — /x^rit = —0.826 into 



For particles with velocity v = c. 
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Fig. 5. Monte-Carlo Simulation of the density profile at an 
oblique shock with compression ratio r = 4, maximal scattering 
angle Af^max ~ 0.1, shock speed Us = 0.4 and inclination angle 
$ = 45° (compare Fig. U and Fig. m. The solid line shows 
the 'bins' of a position measurement of particles, whereas the 
filled dots show a density measured through the flux through 
a surface with constant distance upstream of the shock. The 
flUed square indicates a measurement downstream of the shock. 



the upstream rest frame gives (— /i„it)* = —0.949. At the 
downstream side of the shock the distribution in the dHT 
frame (solid Une of Fig. ||b) is simply divided into parti- 
cles going upstream (— 1 < /x < 0) and those going down- 
stream (0 < /i < 1). The boundary /Ut = between these 
regions is transformed to /i* = —0.446 in the downstream 
rest frame. 



3. 3. Density profile 

From the pitch-angle distributions one can calculate the 
density using Eqs. (5) and (6), which is the integration 
of the pitch-angle distributions in the dHT frame. From 
Fig. we see that the distribution immediately down- 
stream of the shock is ^ 0.38 (open circles in Fig. Hb), 
compared to < 1 for the upstream distribution (open cir- 
cles in Fig. ^). By comparing the upstream and down- 
stream distributions in the dHT frame, it is obvious that 
the downstream integral is less than upstream, and there- 
fore the density is discontinuous at the shock. 

In planar symmetry, the density distribution of accel- 
erated test particles is a function of the distance from the 
shock front. We use the dHT frame in order to compare the 
densities upstream and downstream directly, as in Sect. |[ 
In this frame the shock is stationary (at a; = 0). We nor- 
malise the distance x perpendicular to the shock as the 
dimensionless variable ^ according to 
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Fig. 6. Ratio of the upstream to the downstream density at 
the shock front vs. the intersection velocity of the shock and 
the magnetic field for infinitesimal pitch-angle scattering with 
ASlniax ~ 0.1, compression ratio r — 4 and three different shock 
speeds. The values of n/n' are taken from flux measurements 
at the shock front, exemplified by the filled circle and square 
at ^ = in Fig. |[ 
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where $ is the upstream inclination angle of the magnetic 
field as defined above, and Ky the parallel diffusion coeffi- 



where Vs is the scattering frequency and we have used 
V = c in the normalisation. Note that Ky , while conve- 
nient for the determination of a length scale, does not 
describe the transport of particles for the anisotropic dis- 
tributions discussed here. The mean squared variation in 
H can easily be calculated from Eqs. (9) and (10). Defining 
SfiQ := cos(Ar2,-„„^), then 

{(Afif) = i(i - <5mo)(i - S^,o^^') - i(i - J^^o)' ■ (15) 

Figure || shows the steady-state density in the dHT frame 
for a shock with velocity Us = 0.4, compression ratio r = 4 
and inclination angle <& ~ 45°. The upstream plasma ve- 
locity in the dHT frame is then Us/cos4> = 0.5657. The 
plot shows the density n(^) as a function of the distance 
to the shock ^. Contributions to this density were limited 
to particles from the power law part of the spectrum, indi- 
cated by the fit in Fig. H. Because the question of spatial 
resolution is of crucial importance to our discussion, we 
use two independent methods to evaluate the density. The 
solid line shows the contents of spatial 'bins', where par- 
ticles simply contribute after every time step to the 'bin' 
at their actual position. These 'bins' are located such that 
the shock lies at the border between two of them. An inde- 
pendent way of measuring the density is to count particles 
which cross a plane at a certain position. The count rate 
is proportional to the flux through this plane; to obtain 
the density, one divides this quantity by the relative veloc- 
ity of the binned particle and the plane. The shock plane 
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Fig. 7. Ratio of the upstream to the downstream density at 
the shock front vs. the intersection velocity of the shock and 
the magnetic field for intermediate pitch-angle scattering with 
Af2max = 1-0, compression ratio r = 4 and three different shock 
speeds. The values of n/n' are taken from fiux measurements 
at the shock front. For values of Us/ cos <1> ^ 0.8, the results are 
essentially the same as those shown in Fig. fel 



and other planes at constant ^ are stationary in the dHT 
frame, so that the particle velocity relative to the planes is 
vfi cos $dHT in the upstream region (see footnote || for the 
relation between <I> and ^ciht), where the particle velocity 
V is multiplied by the cosine of the angle with respect to 
the magnetic field fj, and the cosine of the inclination angle 
of the magnetic field relative to the normal of the plane. 
The three filled dots at ^ = 0, 0.5 and 1 show the density 
in the upstream region, whereas the filled square shows 
the density at the downstream side of the shock at ^ = 0. 
In particular, the two values at ^ = are calculated from 
the integration of Eqs. (5) and (6) over the dHT distribu- 
tions shown in Fig. (open circles). The normalisation is 
in each case such as to give unity far downstream of the 
shock. 

Both methods display a discontinuous density profile 
in Fig. H. However, finite spatial resolution means that the 
binning method systematically underestimates the value 
of the density when approaching the shock front from ei- 
ther side. The second method, on the other hand, is a 
precise measure of the density at ^ = 0. For the exam- 
ple of Fig. H (10^ independent particles) we get n/n' — 
3.641 ± 0.004, where n represents the filled circle, and n' 
represents the filled square at ^ = 0. Figure || shows the 
ratio n/n' for various shock velocities and inclination an- 
gles. (The number of particles lies between 6 • 10^ to 5 • 10^ 
and the statistical error of n/n' is < 5%). For $ = (par- 
allel shock) we get a continuous density (n = n') for all 
shock speeds, because here no refiection can occur (see 
Sect.g). If the velocity of the intersection point of the 
shock and the magnetic field exceeds 0.8 of the particle 
velocity (in the case v = c discussed here), the upstream 
density becomes more than 10 times the downstream den- 
sity at the shock front. This ratio increases very rapidly 




Fig. 8. Monte-Carlo Simulation of the density profile at an 
oblique shock (<1> = 60°) with compression ratio r = 3.5 and 
shock speed Ub — 300 km s~^ , for various regions of the parti- 
cle velocity. The normalisation of each plot is chosen to avoid 
overlap. The maximal scattering angle is Ar2,„ax ~ 0.5. The 
filled dots show the upstream density measured from the flux 
at the shock and normalised to the value of the corresponding 
plot far downstream. 



for higher intersection velocities. Our calculations are per- 
formed for test particles, but in reality such particles may 
exert a substantial pressure, which would be important in 
calculations which include the back-reaction of the par- 
ticles on the fiow, as pointed out by 091. However (as 
discussed in Sect. pl\ ), increasing the magnitude of the 
maximal change in pitch angle (AO„j,^) leads to a strong 
reduction of the highest possible upstream densities. This 
can be seen very clearly from Fig. u[ where the maximal 
change in pitch angle is Arj^j,^ — 1.0, and arises because in 
this case the angular distribution cannot develop a sharply 
peaked structure such as that shown in Fig. m For isotropi- 
sation of the pitch angle after each scattering event, the 
ratio n/n' for r = 4 is always smaller than 2 for all in- 
clination angles which are possible for subluminal oblique 
shocks (ms < cos$). 

Finally, in Fig. @ we show the density profile for an 
oblique shock ($ = 60°) with non-relativistic shock ve- 
locity Us = 300 km s~^. We have simulated the accelera- 
tion of test particles at an unmodified low Mach-number 
shock with compression ratio r = 3.5, typical for shocks 
in the solar wind. The pitch-angle scattering is described 
by Eqs. (9) and (10) with Af},„,, = 0.5. This plot shows 
the density for 4 different velocity regimes of the acceler- 
ated particles, which were injected with velocity p/m = 
2.1 • 10~^, just higher than the intersection velocity of 
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shock and magnetic field (us/cos^ = 2 • 10^'^). The den- 
sity is normahsed such that overlap of the four plots is 
avoided. For each plot in Fig. || the filled dot represents 
the density measured at the shock front on the upstream 
side from the flux at the shock plane (compare Fig. [^ and 
the corresponding text). For the lowest energy particles 
shown in Fig. |g, which would correspond to protons with 
a kinetic energy between 12 and 19 keV, the density up- 
stream of the shock exceeds the density far downstream 
by more than a factor of 3, and the ratio of the density up- 
stream relative to downstream at the shock front is given 
by n/n' = 2.81 ± 0.02 at this 'energy band' (note that the 
error represents only the statistical fluctuation) . It can be 
seen from Fig. g that these ratios depend on the velocity 
of the particles. If the particle velocity exceeds the veloc- 
ity of shock and magnetic field intersection by an order of 
magnitude, the distribution becomes more isotropic and 
therefore the density jump tends to disappear, as discussed 
in Sect. and |. 

4. Discussion 

We have presented an analysis of particle acceleration at 
oblique shocks, with special emphasis on the density of 
accelerated test particles. It was shown analytically that 
in general a density jump can occur at an oblique shock 
front. This was done on the basis of Liouville's theorem 
and by integrating the phase space density, without the a 
priori assumption of an isotropic pitch-angle distribution. 
It turns out that for situations in which the pitch-angle 
distribution is non- isotropic, a density peak at the shock 
front can appear. This peak has a discontinuity at the po- 
sition of the shock, if the adiabatic treatment is used. One 
can describe the resulting precursor of upstream acceler- 
ated particles as due to reflections at the shock front. 

We used Monte-Carlo simulations to calculate the 
spectral index and the pitch-angle distribution. These re- 
sults are in good agreement with semi-analytical calcu- 
lations from K&II89. The corresponding density profile 
shows a pronounced discontinuity for a large range of pa- 
rameters, which was also found by 091 (see Fig. 6 therein), 
and our results are quantitatively in very good agreement 
(r = 5.28, Us = 0.3, Af7„,,, = 0.1, u,/cos$ = 0.6 => 
n/n' = 6.31 ± 0.04; and Us/cos$ = 0.8771 => n/n' = 
393 ± 7) . This discontinuity results from the larger den- 
sity of reflected particles as compared to transmitted ones 
than is found in the isotropic case. Especially in the case of 
infinitesimal pitch-angle scattering, the particles undergo 
repeated reflections (by which they are accelerated) before 
they reach a pitch angle at which they are able to cross 
the shock into the downstream region. 

Allowing for a larger maximum value in the change 
of the pitch angle increases the probability of entering 
the loss cone per scattering event, and therefore cross- 
ing the shock from upstream to downstream. This re- 
duces the density contrast (Fig. |^). The effect on the spec- 



tral index is restricted to a small region in the parameter 
Us/ cos $ (the intersection velocity of magnetic field and 
shock) (Fig. 0). The reason for this steepening is that the 
acceleration due to reflections becomes less effective. For 
large pitch-angle scattering where the pitch angle is ran- 
domised after every scattering event, the pile-up effect is 
almost absent, and the minimal spectral index (the flat- 
test spectrum) which can be reached at oblique shocks 
with non-relativistic shock velocities is s > 3.4 (see also 
N&T95). 

In the case of non-relativistic particles accelerated at 
solar system shocks, a significant density peak can occur 
only for highly oblique magnetic fields for particles whose 
velocity is less than an order of magnitude greater than 
the intersection velocity of shock and magnetic field. The 
occurrence of a density peak of accelerated particles at the 
shock front could in principle be detected in situ by space 
observations in the solar system, so that it is important to 
determine whether effects such as modification of the ve- 
locity profile by the pressure of accelerated particles can 
affect our results. Simulations of this n online ar problem 
have been performed by Ellison et al. ( 1996 ). However, 
they detect no difference in the spectral index between 
runs with infinitesimal and large pitch-angle scattering, 
whereas we predict that a difference should accompany a 
high density peak ahead of the shock front, as discussed 
in Sec. 3.1. Very recently, numerical solutions of the trans- 
port equation for mildly relativistic particles at solar sys- 
tem shocks have found a density peak at the shock front 
(Ruffolo [l999| ). 

A pile-up of electrons ahead of the blast-wave associ- 
ated with a supernova remnant could cause a synchrotron 
flux which decreases downstream of the shock instead of 
increasing due to the compression of the magnetic field 
(091). Whereas this could always happen at a shock which 
propagates in an inhomogeneous medium, such a signature 
could be also produced by a shock moving in a homoge- 
neous interstellar medium with an oblique magnetic field 
as a result of a pile-up of reflected particles. 

The synchrotron emission as a function of frequency 
scales as e(i^) ex ni'~" B"'^^, where a = {s — 3)/2. Be- 
cause of the enhanced synchrotron emission due to the 
increasing magnetic field downstream of the shock, the 
observation of the density peak upstream of an oblique 
shock would be possible only for very high obliquities ob- 
served at high angular resolution. For a typical supernova 
shock velocity of Ug = 7000 km s^ , an upstream incli- 
nation angle of $ <; 88° leads to an upstream density (of 
highly relativistic particles) which exceeds by more than a 
factor 10 the downstream density n/n' ^ 10 (see Fig. Bfor 
infinitesimal, and Fig. [7| for intermediate pitch-angle scat- 
tering). At the same time the very flat spectrum produced 
by these particles (s « 3) leads to only a linear enhance- 
ment of the downstream synchrotron emission due to the 
compressed magnetic fleld, {B/B')"'^^ w 1/4 (for r = 4). 
This means that the synchrotron emission upstream would 
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exceed the emission downstream. We estimate that such 
an effect could be resolved in the radio range and for pa- 
rameters typical of a young SNR such as Tycho if the scat- 
tering frequency Vs, defined in Eq. (14), is smaller than the 
gyro- frequency by a factor of roughly 10^. 

Observations of the Tycho supernova remnant have re- 
vealed a density peak in the vicinity of the blast wave 
(Reynoso et al. 1997 ), which cannot be understood by the 
theory of an unmodified parallel shock moving in a ho- 
mogeneous medium. At least one of these assumptions 
must be dropped in a model of that part of the rem- 
nant which shows this feature. A highly oblique shock is a 
possibility, especially if nonlinear effects in the precursor 
cause preferential alignment of the field in the plane of 
the shock, as suggested by recent simulations (Bell, pri- 
vate comm.). However, the observed spectral indices of 
young shell-like supernova remnants are in the range 0.5 ^ 
a ^ 0.8 (Reynolds & Ellison 1992), which corresponds to 
4 ^ s ^ 4.6 and deviates from the standard result for 
compression r = 4 at a highly oblique shock. Spectral in- 
dices in this regime can be produced by stochastic perpen- 
dicular magnetic fields (Kirk et al. 1996^ Gieseler 1998) 
or strongly modified oblique shocks (Reynolds & Ellison 
1992). Whereas in the former case a densitv peak ahead 
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of the I shock does not arise (Gieseler 1998). in the latter 
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case the pile-up of accelerated particles may be an impor- 
tant effect. Alternatively, nonlinear hydrodynamic effects 
of the accelerated particles can create an unstable den- 
sity sp ike do wnstream of the shock, as shown by Jun & 
Jones ( 1997 ), which may be responsible for the enhance- 
ment in the synchrotron emission. 

We have shown that highly oblique shocks can produce 
a pronounced density peak due to a pile-up of accelerated 
particles ahead of the shock front. In Monte-Carlo simu- 
lations, this might appear as a discontinuity or even be 
overlooked, depending on the method used and the spa- 
tial resolution achieved. We point out that the nonlinear 
effects of such a pile- up could be significant, so that it is 
important to locate this effect in simulations which incor- 
porate the reaction of the pressure of accelerated particles 
on the plasma dynamics. 
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